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Abstract. In a recent series of papers [Tl[2l|3], a statistical model that accounts for correlations between 
topological and geometrical properties of a two-dimensional shuffled foam has been proposed and compared 
with experimental and numerical data. Here, the various assumptions on which the model is based are 
exposed and justified: the equiprobability hypothesis of the foam configurations is argued. The range of 
correlations between bubbles is discussed, and the mean held approximation that is used in the model is 
detailed. The two self-consistency equations associated with this mean held description can be interpreted 
as the conservation laws of number of sides and bubble curvature, respectively. Finally, the use of a “Grand- 
Canonical” description, in which the foam constitutes a reservoir of sides and curvature, is justihed. 
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1 Introduction 

Foams, granular materials, biological systems, or glasses 
share the common feature of being out-of-equilibrium sys¬ 
tems from the point of view of standard statistical me¬ 
chanics: either these systems possess a large number of 
metastable states in which they are trapped within a usual 
range of temperature (e.g.: foam or granular materials at 
rest, glassy systems); or a flux of energy is continuously 
injected from non-thermal sources and is lost through in¬ 
ternal friction, allowing them to explore their phase space 
(e.g.: sheared or tapped granular systems, biological cells). 
As these systems involve many individual objects (bub¬ 
bles, grains, cells, or atoms), there is a strong motivation 
to treat them with thermodynamic methods. The pioneer¬ 
ing work of Edwards on granular matter nil] has shown 
how the powerful arsenal of statistical physics can be ex¬ 
tended to out-of-equilibrium systems. Because of the per¬ 
fect space-filling property of a foam {i.e. without gaps nor 
overlaps), and the high deformability of its constituents 
(bubbles), its metastable states are described by a set of 
well established rules |5|[7]. For that reason, foams are very 
good candidates of out-of-equilibrium model systems. 

Recently, a new statistical description has been pro¬ 
posed to describe the structural properties of a two- 
dimensional (2D) dry foam [T]. The model succeeded to 
reproduce, without any adjusting parameter, the correla¬ 
tion between geometry (distribution of bubble size) and 
topology (distribution of side number) for foams with 
moderate size dispersities P1I5] . No semi-empirical laws 
(such as Aboav-Weaire or Lewis laws [5117]) are used in 
this approach; instead, laws that correlate geometry and 


topology properties for individual bubbles and for the en¬ 
tire foam are deduced from the model. Interplay between 
geometrical and topological features of a foam is crucial 
in determining e.g. rheological properties or coarsening 
rate. Although the primary goal of this approach is to 
describe the statistical properties of foams which are uni¬ 
formly shuffled, it may also describe disorders in foams at 
rest, if the foam production process involves uniform shuf¬ 
fling. Furthermore, because of the resemblance of many 
epithelial tissues with densely packed soap bubbles m, 
this approach could be extended to biological tissues. Al¬ 
though first comparisons with numerical and experimental 
data show a remarkable agreement [2115] . this theoretical 
model is based on successive assumptions which need to 
be tested. 

The aim of the present paper is to establish solid foun¬ 
dations to this theoretical framework, by discussing and 
testing the validity of the different underlying hypothe¬ 
ses. In particular, the mean held approach used in the 
model is detailed, and one shows that the associated self- 
consistency equations can be interpreted as the conserva¬ 
tion equations of two scalar quantities. The validity of the 
hypotheses are tested with simulations of shuffled foams. 
Simulations are well suited because they allow to have 
precise measures on the pressure inside every bubble and 
curvature of their sides, which cannot be easily done in 
experiments. 

The structure of the paper is as follows: Sectionj^sum- 
marizes the main results of the model and its success to 
reproduce experimental and numerical data. Section [pre¬ 
views the geometrical, topological and mechanical laws 
that govern the structure of a foam at rest. Because of 
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these space-filling constraints, bubble shapes are corre¬ 
lated. In Section]^ mechanical and thermal shufflings are 
defined and compared, and the postulate of equiprobabil- 
ity of the visited states is discussed. Section presents 
a mean field theory in which correlations between neigh¬ 
bouring bubbles are disregarded, but correlation between 
the geometry (size) of a bubble and its topology (number 
of sides) is taken into account. Within the frame of this ap¬ 
proximation, two invariants are identified in a 2D shuffled 
foam, i. e. quantities that are exchanged by the bubbles but 
whose sum is preserved through the elementary topologi¬ 
cal process of neighbour swapping (T1 event). The mean 
field approximation is then tested with numerical simula¬ 
tions. Section 1^ discusses the validity of the assumptions 
made to use a grand-canonical description: the short cor¬ 
relation length between the bubbles, and the extensivity 
of the state variables. Finally, Sectiondiscusses possible 
extensions of the model to three-dimensional foams and 
other cellular systems. 


2 Main results of the model 

In this section, one briefly recalls the main predictions 
of the model, and its success to reproduce, without any 
adjusting parameter, the correlations between topological 
and geometrical properties of a foam. The model concerns 
soft (or liquid), two-dimensional foams in the low density 
limit, so bubbles have polygonal shapes; their structures 
are governed by space-filling constraints that are detailed 
in Sect. Here density is far below the jamming transi¬ 
tion point HZ] and is not a control parameter of the sys¬ 
tem. One focuses on time scales much shorter than those 
typical of bubble coarsening and coalescence, so the pre¬ 
serve their integrity and size. Strictly, the area Ai and 
pressure Ui of any bubble i depend on the specific ar¬ 
rangement of all the bubbles, and so their values fluctuate 
in a shuffled foam; only the amount of gas that it con¬ 
tains is conserved. However, as discussed in Appendix A, 
the area and pressure fluctuations with foam configura¬ 
tions are usually very small [Tj and so Ai can be equated 
with its average value: Ai ~ A^, where the upper bar de¬ 
notes a time average or, assuming ergodicity, an averaging 
over the foam configurations (ensemble average). Hence, 
each bubble can be unambiguously identified with its area 
rather than the amount of gas it contains. To simplify the 
expressions hereafter, it is convenient to introduce the ef¬ 
fective radius Ri = \JAijir, i.e. the radius of the circle 
which has the same surface area. 

There are experimental and numerical evidence m 
that the geometrical (size) and topological (number of 
sides) features of a bubble are correlated: in a 2D monodis- 
perse foam (one bubble size), most of the bubbles are 6- 
sided. Conversely, in a polydisperse foam, the larger bub¬ 
bles have statistically more sides than the smaller ones. 
Correlations between geometrical and topological quanti¬ 
ties of a foam are characterized by the conditional prob¬ 
ability P{n\R) for a bubble to have n sides, given its size 
R (conditional side number distribution). The model pre¬ 


dicts [T|: 


P{n\R) 
where ^{R) = 


exp 


- &) + ^.n , 
-|^n(n - 6) -I- 


( 1 ) 
is the 


partition function of the bubble, is the minimal num¬ 
ber of sides of a bubble, often taken as equal to 3 in the 
literature, but in principle n^nin = 2 is not forbidden by the 
space-filling constraints. There is no free parameter in the 
model: e ~ 3.72, while (3 and /x are implicitly related to the 
(known) size distribution P(R) through: d\x\S/df3 = 0 
and i9 In S '/dgL = 6Af, where S is the partition function of 
the entire foam, defined through 


( 2 ) 


No assumption is made on the shape of the size distribu¬ 
tion P{R)', it can be unimodal or not. For a given P{R), 
the equations of the model can be solved numerically. The 
model predicts an order-disorder transition at a small, but 
finite, size dispersity |3]. Below this transition point (crys¬ 
tallisation threshold), all bubbles are 6-sided cells: the 
foam is topologically ordered. For a bidisperse foam, the 
transition occurs at the large-to-small size ratio 7/5 = 1.4. 
Above this point, analytical approximations can be de¬ 
rived: the truncated and discrete normal distribution (IH 
is evaluated by treating n as a continuous variable jHI^ • 
This approximation yields analytic expressions for /? /x: 


((«-■) 

(3) 

-i=3((i?)(i?-i)-l), 

(4) 


Here the brackets denote an averaging over the bubbles 
that compose the foam: for any quantity Xi attached to 
every bubble i, (x) — ■ Introducing p{R,n), 

the proportion of bubbles with size R and number of 
sides n, this average can also be expressed as (x) = 
n} n)x{R, n). It should not be confused with the 
time (or ensemble) average defined above. 

Eqs. ^ -0 allow to express the average side number 
of a bubble with size R: 

The topological disorder anl{n) = ^ (n^) — {n)"^ j(n) can 
also be expressed as a function of the (known) moments 
of the size distribution: 

(^) + ( 6 ) 

In principle, the analytic expressions (§,@0,0 are 
valid only for sufficiently high poly dispersity of size. In 
practice, these are very good approximations down to the 
crystallisation threshold [3|. 

Predictions on individual [Eq. 01 and global [Eq. 
01 topology-geometry correlations agree well with ex¬ 
perimental and numerical data of foams that are shuf¬ 
fled, either thermally or mechanically m\- Nevertheless, 
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because of the very different nature of these two shuf¬ 
fling mechanisms, one does not expect to observe such 
agreement for higher shear rate or higher temperature, 
nor for all statistical properties of the foam. In particular, 
one expects the two shuffling mechanisms present different 
barrier-crossing statistics mM- 


3 Space-filling constraints in a foam at rest 

A 2D foam is a partition of the plane without any gaps or 
overlaps, and its structure must obey a set of well estab¬ 
lished physical, geometrical, and topological constraints 
|F1I7|. The physical constraints follow from the mechanical 
equilibrium throughout the system: the balance of film 
tensions at every vertex yields Plateau’s laws: the edges, 
or Plateau borders, meet in threes at a vertex, at an an¬ 
gle of 120° with each other. The balance of normal forces 
across each film yields Young-Laplace’s law: every edge 
is an arc of circle, whose algebraic curvature Hij is pro¬ 
portional to the pressure difference between the adjacent 
bubbles i and j: 

7T, - Ui .. 

Ktj=-Kj^ = - - -, (7) 

where 7 is the 2D film tension, and 11 i and 11 j are the 2D 
pressures in bubble i and j, respectively (by convention, 
Kij > 0 when the center of curvature is outside the cell 
i, i.e.'. when 11 j > Ui). It must be emphasized that at 
microscopic scale the interface is not smooth, because of 
thermal fluctuations. In fact, Young-Laplace’s law relates 
the mean curvature to the mean bubble pressures. 

There is also a constraint on every bubble because of 
the fixed amount of gas that it contains. This constraint 
is often replaced, for simplicity (see Appendix A), by a 
constraint on every bubble area (incompressible gas). 

Euler’s rule is a topological constraint which relates 
the number of bubbles M, with those of Plateau borders 
Npb, and vertices Ny in a foam that covers all the 2D 
space: 

M-Npb + Ny = X, (8) 

where x is the Euler-Poincare characteristic {e.g. x = 0 
for a toroidal space, x = 2 for a spherical one). 

Finally, the Gauss-Bonnet theorem relates geometry 
and topology of every bubble i: 

where iTl(j) represents the set of the Ui bubbles that neigh¬ 
bor bubble i, and is the length of the edge common to 
bubbles i and j. 

The constraints listed above substitute for the steric 
repulsion (excluded volume) that lies in hard-core granu¬ 
lar materials: as a consequence, there is no limit on the 
number of neighbours that a bubble can have (while in 
a monodisperse disk packing, no more than 6 neighbours 


can fit around a disk because of steric exclusion). Curva¬ 
ture and length of the sides can adapt to satisfy Young- 
Laplace’s and Plateau’s laws; this induces an increase of 
pressure difference between the central bubble and its 
neighbours. Geometry and topology of the bubbles are 
then correlated through these space-filling constraints. 


4 Shuffled foams 

4.1 Mechanical and thermal shuffling of foams 

From the point of view of standard thermodynamics, a 
liquid foam at rest is an out-of-equilibrium system: it has 
many metastable states, but the energy required to jump 
from one metastable state to another is much higher than 
the thermal energy kpT within a usual range of temper¬ 
ature. Thus, in absence of coarsening and coalescence, it 
takes an indefinitely long time to escape from an energetic 
valley. Actually, the temperature required to trigger a T1 
event is far above the vaporizing temperature of the liquid 
phase, which makes thermal shuffling inaccessible exper¬ 
imentally. Numerically, on the other hand, this mode of 
shuffling is easy to handle. This can be done with Monte 
Garlo simulations (see Appendix B). 

Mechanical shuffling also confers dynamics to the foam 
and can be achieved both numerically and experimentally. 
However, such a system is still out of thermodynamic 
equilibrium: energy is injected (in the form of mechani¬ 
cal work) to the system by a non-thermal source, which 
is then dissipated during the rearrangements associated 
with T1 events |21| . 


4.2 Gentle shuffling of foams 

Foams can explore their phase space by various shuffling 
mechanisms. One restricts this study to gentle shufflings, 
defined as follows. 

For mechanically shuffled foams, one assumes that 
shuffling is very slow when compared to the typical re¬ 
laxation time associated with a T1 event: the transient 
states can then be neglected in the ensemble of accessible 
states: foam passes instantaneously from one metastable 
state to another. Moreover, one considers only shuffling 
mechanisms where the bubble areas, and hence the total 
foam area, are preserved. Because of the high deformabil- 
ity of its constituents, such a restriction is always possible 
for a foam, unlike for jammed hard granular systems | 22 | . 
One also supposes that T1 events occur homogeneously 
within the foam. This assumption certainly restrains the 
mode of shuffling: near the box boundaries (shear banding 
phenomenon) |241I25] . Gonversely, a random, non periodic, 
shuffling should satisfy this hypothesis. It seems that such 
an experiment has never been done yet. 

Thus, for both shuffling mechanisms, the structure of 
the foam evolves though T 1 events exclusively, and bubble 
preserve their integrity. 
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4.3 Equiprobability hypothesis 

By analogy with the fundamental postulate of statisti¬ 
cal mechanics |26] . one assumes that, under a gentle and 
homogeneous shuffling, all the accessible microstates of a 
large foam {Af ^ 1) are visited with equal probability. 

One may wonder whether or not the equiprobability 
hypothesis must be restricted to the states with same en¬ 
ergy, as it is postulated for thermal systems at equilibrium. 
Actually, the answer depends on the origin of the shuffling. 

For foams which are mechanically shuffled, the answer 
to this question is still a large field of investigation: as 
previously mentioned, such system is out-of-equilibrium 
from the point of view of conventional thermodynamics 

miMm- 

Fig. [3 shows the fluctuations observed for different 
quantities attached to a single bubble i within a thermally 
shuffled foam, using Cellular Potts Model simulations (see 
Appendix B). As expected (see Appendix A), variations of 
the bubble area are much smaller than those of the pres¬ 
sure difference |7Tj - where 77^^) = ^j/rii 

is the instantaneous average pressure in the bubbles that 
belong to Df(*), the neighbourhood of bubble i. 

Energy fluctuations are also very small. Almost no en¬ 
ergy is exchanged between bubbles, which justifies that 
the Boltzmann weighting can be neglected in a gentle 
thermal shuffling. This is also consistent with previous 
studies that show that the difference of energy between 
metastable states are less than 2% for a 2D foam PTII^. 
In the mean-field approximation detailed in Sect. bub¬ 
bles have constant perimeters: the energy of every bubble, 
and hence of the entire foam, are fixed. 

For these reasons, one can safely disregard the energy 
weighting of the probability distribution, and not distin¬ 
guish thermally from mechanically shuffled foams. Indeed, 
in spite of the profound difference of these two modes of 
shuffling, experimental and numerical data collapse ex¬ 
tremely well mg in the gentle shuffling regime. 

4.4 Degrees of freedom of a gently shuffled foam 

To perform a statistical description of the foam structural 
properties, one needs to enumerate the microstates which 
correspond to the same macroscopic state. A detailed de¬ 
scription giving positions and momenta at the atomic scale 
is too cumbersome to handle, and one must use instead a 
coarse-grained description in which the state of the foam 
refers to its structure (topology and geometry) at the bub¬ 
ble scale. For instance, the complete knowledge of the con¬ 
tour of every bubble fully describes a metastable state • 
Actually, contours are not independent in a foam at rest, 
because of the space-filling constraints. The minimal num¬ 
ber of degrees of freedom that are required to describe a 
foam metastable state is still an open question |28l[Ml[Ml 
mm- Probably the most comprehensive study on this 
domain has been carried out by C. Moukarzel which has 
identified a set of 4M independent variables to fully de¬ 
scribe any metastable configuration of a foam that tiles 
the plane [3g. However, these variables have no direct 
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time [MGS] 

Fig. 1. Typical fluctuations for different quantities attached 
to bubble i: area Ai (purple curve), energy £i (green curve), 
and pressure difference \IIi — (blue curve), where 77^^) = 

X)jg(n(i) 77j/ni. Curves are shown on two different graphs be¬ 
cause of the difference in amplitude of the fluctuations. Bubble 
perimeter is not shown, as it is barely distinguishable from the 
energy curve. Topological charge qt{ni) = 6 — ri; (red curve) 
is also drawn to indicate the T1 events in which bubble i is 
involved. 


physical interpretation. Here instead, one assumes that 
a microstate L of the foam is correctly described by the 
number of sides ni, location and pressure 11 i of ev¬ 
ery bubble: L = i = 1, 2,..., Af}. Although 

pressure fluctuations are small, the fluctuation of pressure 
difference between adjacent bubbles can be large (see Ap¬ 
pendix A). That is why they must be taken into account 
in the description of a microstate. 

, the microstate of the foam L is given by the mi¬ 
crostate ii = {ri,ni,IIi) of every bubble i, which is con¬ 
venient to apply a mean field description. However, this 
choice is not crucial for the rest of the paper: the statistical 
description presented here only requires that a microstate 
of the foam is adequately described by the variables {ui}, 
plus possible other degrees of freedom which must be in¬ 
dependent of the {ui}. 

When a bubble i is directly involved in a Tl, its asso¬ 
ciated variables (ri,ni, Ili) vary together, inducing corre¬ 
lations between them. However, bubble positions {r^} are 
continuous variables and have some uncertainty attached 
to them jSl]. Thus, after few Tl events, and (ni,7Ti) 
can be fairly treated as independent variables: at a given 
value of {ni^Ui) correspond many different values of r^, 
and vice-versa. Note that the same kind of assumption is 
done in the statistical mechanics of a gas (Stosszahlansatz 
or molecular chaos hypothesis m) to justify the inde¬ 
pendence of positions and momenta of the particles. In 
addition, thanks to the high deformability of the bubbles, 
(jii^Ui) are almost independent of the relative positions 
of the N — 1 other bubbles {r^ — r^, j ^ i}. Therefore, 
the variables {r^} and {ni,i7i} can be treated as two sets 
of variables that are independent from each other. 
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It is shown in the next Section that, within a mean 
field description, variables {r^} and {rij} are sufficient to 
properly describe a foam microstate. 

5 Mean Field Approximation 

Because of the space-filling constraints (see Sect.j^, bub¬ 
bles are correlated objects : a change in the geometry or 
topology of a bubble affects the geometry and topology 
of the other bubbles. Therefore, the number of accessi¬ 
ble states cannot be simply obtained by enumerating the 
possible distributions of the N sides over the JV bubbles: 
many of these distributions are not compatible with the 
space-filling constraints. A perfect description taking all 
these constraints into account is illusive. Instead, a mean 
field approach is proposed, in which the fluctuations of the 
neighbourhood of a bubble with given size R and given 
number of sides n are neglected. To legitimate this ap¬ 
proximation, one shows below that in most cases these 
correlations are in fact short-ranged. 

5.1 Screening of bubble correlations 

One discusses how the bubble microstates £i = (rij, r^, 77^) 
are correlated to each other. One should first study how 
they are altered by a unique T1 event in the foam: clearly, 
only the four bubbles involved in the T1 have their rii that 
vary. The pressure field is a priori modified on a longer 
range: previous studies EH1I35] have shown that the ex¬ 
tra pressure due to a topological charge qt{n) (a bubble 
with n 7 ^ 6 sides) within a regular foam of hexagonal bub¬ 
bles decreases logarithmically (as for the electric potential 
generated by a point charge in a 2D dielectric medium). 
However, a large foam must contain in same proportion 
positive {qt = 6 — n > 0) and negative (g* < 0) topological 
charges. For instance, a T1 event in a hexagonal foam pro¬ 
duces two positive and two negative topological charges 
simultaneously. Their antagonist effects on the pressure 
field screen their correlations on a typical lengthscale of 
one bubble size. This “screening” effect is further enhanced 
in a disordered foam |3S]. 

Displacements are also very localized: Figure shows 
the superposition of the pictures of a dry foam before and 
after a T1 event, either experimentally or numerically. Es¬ 
sentially the four bubbles directly involved in the T1 event 
have their positions and geometries modified; the rest of 
the foam structure is unperturbed. Thus, as long as the 
bubble i is not involved in a T1 event, it remains almost 
in the same microstate ip. its number of sides rii is fixed, 
and the change of its position and pressure Ui are im¬ 
perceptible. Precisely, the fluctuations of position Svi due 
to Tls elsewhere in the rest of the foam are small There¬ 
fore, as long as T1 events remain diluted, correlations be¬ 
tween bubble microstates £i are short-ranged: the variables 
Ui) of a bubble are primarily correlated to those 
of its first neighbours. It can be noticed that correlations 
of bubble displacements are also short-ranged. This sit¬ 
uation is very different from what is observed in dense 



Fig. 2. Superposition of the metastable states of a foam before 
and after a T1 event; (a) experimental foam, excerpt from [36]; 
(b) numerical foam (Surface Evolver), excerpt from |35| . 


packings of hardcore particles: because of steric exclusion, 
the positions of the particles are not independent, yield¬ 
ing macroscopic correlation length of the displacement. 
Direct measures of bubble displacements in partially wet 
foams (liquid fraction ~ 10%) confirm the short corre¬ 
lation length of bubble displacements, typically of a few 
bubble radii IS71I3SI, and hence much lower than what is 
observed in dense granular packings of hardcore particles 

m- 

5.2 Isotropic (regular) bubble model 

The proposed Mean Field Approximation (MFA) consists 
in neglecting the fluctuations of the neighbourhood of each 
bubble with given size R and number of sides n. In other 
terms, the geometry of such a bubble is equated with the 
one obtained from averaging over all the foam configura¬ 
tions for which the number of sides of this bubble is n. 
Since all the constraints listed in Sect. |3]are satisfied for 
every configuration L of the foam, one looks for an aver¬ 
age bubble geometry that satisfies these constraints too. 
Assuming that the foam is homogeneous and isotropic, 
the bubble is surrounded by a uniform foam (identical 
bubbles): all sides are equivalent and the geometry must 
correspond to the geometry of a regular n-sided cell, for 
which sides are identical in shape and length (see Fig.[^. 

The geometry of such regular cell is entirely deter¬ 
mined by its size and number of sides EHiini. In par¬ 
ticular, the curvature of each side is ^/Tr{n — 6)/ (3e(n)i?) 
(by convention, the curvature is chosen positive when its 
center is outside the bubble). e(n) is the elongation of 
the cell (ratio of perimeter to square-root of area). It is a 
very slowly decreasing function, lying between e(2) ~ 3.78 
and e(oo) ~ 3.71. In fact, 2-sided cells are seldom ob¬ 
served in gently shuffled foams, and e(3) already drops 
down to ~ 3.74. Exact expression of e(n) is given in refs 
[T1I281I13| . Hereafter, this weak dependence is neglected, 
and one notes e = e(6) ~ 3.72. It can be noticed that 
this approximation implies that the perimeter of a bub¬ 
ble does not depend on its number of sides; hence, the 
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(a) (b) 


Fig. 3. (a): Real two-dimensional foam in one of its metastable 
states, (b): Mean-held approximation: each bubble i is a regular 
cell surrounded by a uniform and Isotropic foam, resulting form 
the averaging of the positions and side numbers of the N — 1 
other bubbles. 


energy of every bubble, and thus of the entire foam, is 
fixed. The restriction of the equiprobability hypothesis to 
the configurations of same energy discussed in Sect. |4.3| is 
superfluous within this MFA. 

The regular bubble (rb) model presented here assumes 
that the difference of pressure between a bubble i and its 
neighbourhood depends primarily on the size and number 
of sides of that bubble, and not on the size polydispersity 
of the foam. This may fail for very large size-dispersity 
values. 

It will be useful to introduce also the bubble curvature 
of a n—sided bubble, defined as the sum of the algebraic 
curvatures of its sides. For the regular bubble, its expres¬ 
sion is: 

/u V^n(n-6) 

The side curvature of such a bubble is just Krb{R,n)/n. 


5.3 Self-consistency equations 

By construction, a regular bubble satisfies most of the 
constraints listed in Sect.[^ Only two sorts of constraints 
are not automatically satisfied in an idealized foam where 
all the bubbles have regular shapes: a global constraint 
caused by the Euler’s formula [Eq. ®], and local con¬ 
straints caused by the Young-Laplac^ law [Eq. 01 at 
every Plateau border. In the following. Consider first the 
rather academic case of a foam with no boundary, i.e. 
that covers entirely the 2D space. Since all the vertices of 
such a foam are 3-connected, Eq. (|^ yields the first global 
constraint: 

AT 

= 6 (A/'-x). (11) 

i=l 

To derive the second global constraint, it can be first no¬ 
ticed that since a bubble with rii sides belongs to the 
neighbourhood of rii different bubbles in an unbounded 
foam, pressures in adjacent bubbles must obey the follow¬ 
ing identity: 

(12) 


Equations (11 1 and (12) are satisfied for every configu¬ 
ration L of the foam. Therefore they must also be true on 
average. Let pl he the probability for the foam to be in 
configuration L, and Pi{n) =, the probability that bubble 
i has n sides 

(13) 


Averaging of Eqs (111 and (12 1 over all foam configu¬ 
rations {L} leads to: 


Af 

^ ^ p^{n)n = N, 


i=l n>n„i„ 

= 0 , 


(14) 

(15) 


where N = 6 (A/' — x) is the sum of side numbers. On the 
other hand, the MFA presented above consists in neglect¬ 
ing the fluctuations of the neighbourhood for every bubble 
i of given number of sides np. 


Neglecting these fluctuations in Eq. (12 1 yields: 


(16) 

(17) 


for every foam configuration. Introducing P{n\R = Ri), 
the proportion of bubbles with size Ri that have n sides, 
this equation rewrites: 


Af 

E E P{n\R = R,)nAn,{{Ru},n) = h. (18) 

As explained in Sect. |5.H correlation lengths are short 
(compared with the system size) for a foam far from the 
crystallization transition. By design, correlation length is 
zero in the MFA. Averaging over all the configurations of 
bubble i (time or ensemble average) can then be equated 
with averaging over all the bubbles that have same size 
than bubble i within the foam (spatial average), i.e.: 
Pi{n) = P{n\R = Ri). Then, Eqs (15) and ( |I^ are iden¬ 
tical. 

This result shows that the MFA consists in replacing 
the local constraints [Eq. 01 by a unique global constraint 
[Eq. (18)[. This reduction of the number of constraints is 
accompanied by a reduction of the number of degrees of 
freedom necessary to describe a foam microstate L, as 
the regular bubble model relates pressures {Pi} and side 
numbers {rii}. 

The mean bubble curvature is defined as 


(19) 

Using the geometry of a regular bubble [Eq. @ 1 , 
Ki({i?fc},ni) = Ki.b(i?i, rii), and Eq. (15) becomes: 


AT 


Pi{n)K^b{Ri,n) = h. 


( 20 ) 


i—1 n>r 


For a foam, 7 is a constant that has been dropped out. 
However, for the extension of the model to other cellular 
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systems, like epithelial tissues, interfacial tension can have 
different values m, and the regular bubble model must 
be adapted accordingly. 


Equalities (141 and (201 hold for foams with no bound¬ 
aries. They can be easily extended to the case of a free 
bubble cluster (he. a foam that does not cover all the 2D 
space), or a set of individual bubbles within a larger foam. 
Assuming again Pi{n) = P{n\R = Ri), one as: 

N 


Y, Y Pin\R = R^)n = N, ( 21 ) 

i=l n>n„i„ 

Y Y P{n\R^R.)=K. ( 22 ) 


i — 1 n>nmin 


Eqs ( [2T| ) and ( 22 ) constitute self-consistency equations of 
the MEA. A foam with no boundary can be seen as the 
special case where N = 6{Af — x) and K = 0. Otherwise, 
N and K are fluctuating variables. . Clearly, N and K are 
independent variables: there are many ways of distributing 
N sides to the Af bubbles, leading to different values of K 
(and vice-versa when Af ^ 1). 

For a free cluster, adjacency of bubbles imposes that 
sides cancel in pairs, except for those at the boundary of 
the cluster. The curvature K is then reduced to the sum 
of the curvatures of the sides that belong to the outer 
boundary B: 


K = = (23) 

iGB zGB 


where TTq is the surrounding pressure. Values of K and 
N (and therefore of K and N) are restricted within a 
certain range: K < 0 since the pressure in any bubble 
of the cluster is higher than the external pressure |32| . 
Combining Euler’s rule [Eq. ([^] and Plateau’s law gives 
Npb = 3{Af 1 — x) for a cluster. Since two adjacent bub¬ 
bles share the same Plateau border, each Plateau border 
in the bulk contributes to two bubble sides, while each 
Plateau border at the outer boundary contributes to only 
one bubble side. Thus, N < 6{Af -I-1 — x)- 

Self-consistency equations and are also valid 
for a foam that fills a container; in that case N fluctuates, 
while K = K = 0. 


Fluctuations of N and K around their respective aver¬ 
age values are characterized by SN/N and SK/K, where 

SX = \J is the standard deviation of the vari¬ 
able X, calculated over all foam configurations ((5A must 
not be confused with the standard deviation ux calcu¬ 
lated on bubble population). It is shown in Appendix A 
that these fluctuations become negligible as the number of 
bubbles Af increases. However, 6NjN and 6K/K scale as 
Af~°‘ with exponents a > 0 that may differ from the value 
1/2 characteristic of additive variables pS]: indeed, a de¬ 
pends on the boundary conditions considered (free clus¬ 
ter, unbounded foam, or set of separate bubbles). Since 
fluctuations vanish in the limit Af —>■ oo, the postulate 
of equiprobability presented in Sect. |4. 3] can be rephrased: 


all the accessible states with same values of total curvature 
(= K) and number of sides (= N) are equally probable. 


5.4 Interpretation: invariants and conservation laws 


Self-consistency equations pT| and (221 can be inter¬ 
preted as the conservation laws of the numbers of sides Ui 
(or of the topological charge qt{ni)) and the cell curvatures 
KYh{Ri,ni), respectively. Since only the four bubbles di¬ 
rectly involved in a T1 event have their and Krb(.Ri, ni) 
which vary, these quantities are conserved locally. A tight 
parallel can then be drawn between the dynamics of a 
shuffled foam and of an ideal gas, as illustrated in Fig. 
1^ T1 events in a shuffled foam play a role equivalent to 
the collisions between molecules in a gas. As for the mo¬ 
mentum and energy of point particles, the total number of 
sides and bubble curvatures of the four bubbles involved in 
a T1 event are preserved, i.e. they are the invariants asso¬ 
ciated with the dynamics of the foam. For a real foam, the 
space-filling constraints listed in Sect. [^generally imply a 
slight change of geometry of the surrounding bubbles as 
well. For instance, a T1 swap event in an initially regu¬ 
lar hexagonal foam must induce a change of curvatures of 
sides in the neighbouring bubbles in order to satisfy the 
120° angle condition). 


collision (ideal gas) 



T1 event (2D foam) 



Fig. 4. analogy between a collision in an ideal gas and a T1 
event in a 2D shuffled foam, playing the role of elementary 
“thermalization” process. 


5.5 Validity of Mean field approximation 

Mean field description requires that shuffling induces ho¬ 
mogenization of the foam (no segregation of bubbles with 
same size), such that the number of sides of a bubble is 
decorrelated from its position within the foam. Experi¬ 
ments |42] and simulations |43] support this hypothesis. 

Mean field theories are known to work better when the 
number of neighbours increases, as their fluctuations av¬ 
erage out. Therefore, one expects the model to give better 
results with bubbles having many sides {i.e. with larger 
bubbles, statistically). More generally, the mean field ap¬ 
proximation consists in neglecting the fluctuations of bub¬ 
ble curvature. Thus, the approximation is expected to fail 
near the crystallization transition predicted by the model 
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j3]. If this transition does exist, one expects to see an 
increase of T1 events close to it. Although this order- 
disorder transition is well identified in hard granular mate¬ 
rials |4ll4^l46ll47] . it seems that it has never been studied 
numerically or experimentally in foams. 


In summary, the mean field theory consists in ne¬ 
glecting the fluctuations of neighbourhood: ~ 

Uj) depends essentially on rii and Ri, and barely on the 
specific foam configuration; and ii) assimilating the mean 
bubble geometry with the geometry of a regular bubble. 
These assumptions are tested with Potts simulations (see 
Appendix B): Fig.l^shows the values of RiKi{{Ri.}, n) in 
a tridisperse foam for different values of n. Tridispersity 
is a good compromise to have several bubble sizes while 
keeping a large number of bubbles of same size. To avoid 
superposition of hundreds of curves, RiKi{{Rk},n) is cal¬ 
culated over all the bubbles having same size Ri. For a 
given n, the values corresponding to the three different 
bubble sizes collapse very well, what confirms that this 
quantity depends on n but not on Ri. Moreover, its vari¬ 
ation with n is in very good agreement with the regular 
bubble model [Eq. @1- In the inset is plotted the dimen¬ 
sionless side curvature Riki/n for the same data: they are 
aligned on a straight line passing through 0 for n = 6, 
in agreement with the regular bubble model [Eq. (lOl]. 
It can also be noticed that the relative amplitude of the 
fluctuations decreases as n increases, as expected. 



Fig. 5. Dimensionless mean bubble curvature Rk vs. 
n in a tridisperse foam. Orange symbols: small bubbles 
(area= 200 pixels^)-, green symbols: medium bubbles (area= 
400pixels^); blue symbols: large bubbles (area= 600 pixels^). 
Proportion of small, medium and large bubbles are 0.319, 0.335 
and 0.346, respectively. Red curve: regular bubble model [Eq. 
(101[. Error bars and shaded area indicate the fluctuations 
(standard deviation) of Ri around its (par¬ 

tially) averaged value value Riki{{Rk}, n). Inset: dimensionless 
mean side curvature Rk{R,n)/n for the same data. 


Figure l^shows the same quantities compiled for foams 
with various size dispersities (bidisperse, tridisperse, and 
normal distributions of bubble areas). Once again, it shows 
a ve ry good agreement with the regular bubble model [Eq. 
(10l[. In particular, for a given value of n, the data cor¬ 
responding to different foam polydispersity collapse very 
well, which confirms that the mean bubble curvature does 
not depend on the foam polydispersity, on the range of 
polydispersity tested {a a/{A) < 0.8). A linear fit (see 
inset) a(n — 6) to the mean side curvature yields a = 
0.157 ± 0.004, which is very close from the value obtained 
in the regular bubble model (a = y/n/{3e) ~ 0.159). 



Fig. 6. Dimensionless mean bubble curvature Rk vs. n for 
foams with various size dispersities: each point corresponds to 
an averaging over the bubbles with same R and same n within 
a given foam. Each foam contributes several points. Error bars 
represent the Standard Error of the Mean (SEM). Red curve: 
theoretical model [Eq. ( |10[ )). Inset: dimensionless mean side 
curvature Riki/n vs. n for the same data. Red line: linear fit 
y = a{n — 6), with a = 0.157 ± 0.004. 


It must be emphasized that although the MFA neglects 
the anti-correlations between the number of sides of neigh¬ 
bouring bubbles - often modelized by the Aboav-Weaire 
law [Hl F7ll39Ll4nil¥r| - it takes into account the correlations 
between the size and number of sides of every bubble. 


6 Ensemble equivalence 

6.1 Position of the problem 


In the mean held approximation, a microstate of the foam 
is specified by the position and number of sides of every 
bubble {L = {ri,ni}). Since the {ri} an d {m } can be 
treated as independent variables (see Sect. 4.4), the num¬ 
ber of accessible states factorizes as 17 = 1T^17(‘), where 
17(9) is the number of configurations of the geometrical 
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variables {ri}, while 17^*) is the number of arrangements 
of the topological varia bles {rii} tha t are compatible with 
the global constraints (21) and (22). is then a func¬ 
tion of the two state variables K and N. The other state 
variables are Afi,... ,Afp, the number of bubbles with re¬ 
spective sizes i?i,..., Rp, (or equivalently, Af and the bub¬ 
ble size distribution P{R)). For the sake of clarity, only K 
and N are made explicit hereafter. 


In order to study the correlations between bubble size 
and side number, only 1?^*^ needs to be evaluated. Let 
[Ky, N_y\K, N) be the number of microstates of the 
foam for which the curvature and number of sides of a 
subset of bubbles are and N^, respectively. The 
probability that has sides and curvature Ky, given 

their total values N and K in the foam, is then: 


P{K^,N^\K,N) 


{K^,N_s^\K,N) 

{2W(K,N) 


(24) 


However, the evaluation of 17*^*^ is not trivial (its deriva¬ 
tion, above the crystallization, is reported on the Ap¬ 
pendix C) and does not yields a simple expression for 
P(Ky,Ny’\K,N). Instead, a “grand-canonical descrip¬ 
tion” is proposed: the “small” subset of bubbles ^ within 
the foam exchange sides and curvature with the rest of 
the foam through T1 events. The rest of the foam then 
constitutes a reservoir ^ of sides and curvature. However, 
the micro- and grand-canonical descriptions are equiva¬ 
lent for large foams only if the system and the reservoir 
are weakly coupled, that is, if the two following conditions 
are fulfilled: 

i) short-range correlations between the bubbles to ensure 
the statistical independence of =5^ and ^ [iHiiin]. As dis¬ 
cussed in Sect. |5.1[ this is usually true in a real foam, ex¬ 
cept near the cristallization threshold where numerous Tls 
may occur. By definition, these correlations are neglected 
in the MFA. Therefore, 17^*) i^N^, K^\N, K') factorizes as 


{K^,N^\K,N) = {K^,Ny) x 

Q%\Ksi = K-K^,N^ = N-N^), (25) 

where and 17^^ refer to y and when isolated, re¬ 
spectively. 

ii) extensivity of the state variables that define a foam 
macrostate. Therefore, one must check the extensivity of 
N and K. This notion is often (mis) interpreted as a neg¬ 
ligible contribution of the excess (or coupling) term when 
two systems (with at least one of macroscopic size) are put 
together [5IJ1IM] . As an illustration, the energy of Hamil¬ 
tonian systems is extensive when the interactions between 
the constituents are short-ranged. Hence, the thermody¬ 
namic limit can be properly defined. The notion of reser¬ 
voir is properly defined as well, i.e. its effective tempera¬ 
ture and chemical potential are not affected by the pres¬ 
ence of the system in contact. 


6.2 Extensivity and thermodynamic limit 

The scalings of N and K with Af depend whether one con¬ 
siders a set of individual bubbles (taken within a larger 
foam), or a cluster of bubbles. Clearly, in the first case, 
both variables scale linearly with AA. When bubbles are 
clustered, their is a nonlinear contribution of the outer 
boundary to the number of sides. From Euler’s rule [Eq. 
(|^] and Plateau’s law, the number of edges (Plateau bor¬ 
ders) of a cluster is directly related to its number of bub¬ 
bles: Nph = 3(A/'-l-l —x). Since two adjacent bubbles share 
the same Plateau border, the number of sides is twice the 
number of Plateau borders, except for those at the outer 
boundary. Thus, the mean number of sides can be written 
N = dlAA + 1 — x) ~ cVaA, where c is a positive constant 
which depends on the size distribution within the cluster. 

While for N the contribution of the outer boundary 
becomes negligible in comparison with the bulk term as 
J\f increases, this is the one and only contribution for 
K (since side curvatures cancel in pairs in the bulk): 
K = —d'/JA/{R), where once again d is a positive con¬ 
stant which depends on the size dispersion within the clus¬ 
ter. For a foam enclosed in a container, d — Q. For an un¬ 
bounded foam, c = d = 0. For a free cluster, an order of 
magnitude for c and d can be estimated considering a cir¬ 
cular cluster (radius Re) of identical bubbles with radius 
R, as shown on Fig. The perimeter of this cluster can 
be approximated as 271 Re ~ Ni,{2R), and its surface area 
as ttR^ ~ N{77R?), where Aj, is the number of sides at the 
boundary. Therefore, c = Nb/VJA ~ tt. The radius of cur¬ 
vature Kq ^ of a side that belongs to the boundary tends 
to 2R when Af ^ 1. Therefore d = Nb/yOA kqR ~ 7r/2. 



Fig. 7. Boundary of an idealized circular cluster made of iden¬ 
tical bubbles. 


Because of the boundary contribution, N and K are 
not additive quantities: when two shuffled clusters and 
A/i - with respective mean numbers of sides N-y, AA^ and 
respective mean curvature Ar*j*, ^ mixed together, 
their total number of sides A^iu 2 and total curvature A'iu 2 
are not equal to the sum of their values when they were 
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isolated. As an illustration, consider two identical clusters, 


i.e. made with the same number of bubbles {Afi = 
same distribution of bubble size (hence ci = C2, 
and (i?)i = (77)2), Then 

A/2) and 
di = d2 

]V°=]V° = 6(A7i)-ciV(aT 

(26) 

and 

A^iu 2 = 6(2A/i) — ci\/2A/i. 

(27) 

The excess number of sides, defined as A^i2 = Ki\j 2 — K\ — 
N 2 growths sublinearly with Afi : 

A/12 = ci(2 - y2)/A)6. 

(28) 


Therefore, this contribution becomes negligible when Afi 
increases. The mean curvature for each individual cluster 


Kl=Kl = -di.M/{R)i, (29) 

and for the mixed clusters: 

Kiu2 = -di^/Wi/{R)i. (30) 

The excess curvature K 12 = A'iu2 ~ k'i — k\, is here 
of the same order of magnitude of the curvatures of the 
individual clusters: 

(Ai2 = di(2-v^)V((^/(i?)i. (31) 


Nevertheless, for separate bubbles as for clustered bub¬ 
bles, N and K are extensive quantities, meaning that 
N/M and KjM tend to finite limits as Af —t c» |52| . 
Moreover, fluctuations of N and K are always negligible 
for large systems: SN(N —>■ 0, SK/K 0 when Af —t 00 
(see Appendix A) These results together show that the 
thermodynamic limit exists for any large set of bubbles. 
For clustered bubbles in particular, the limits of N/Af 
(~ N/J\f) and K/J\f (~ K/Af) are universal, i.e. they 
do not depend on c and d, or on the specifi c bu bble size 
distribution: N/J\f —>■ 6 and K/Af —>■ 0. Eqs (211 and (221 
can then be rewritten as: 


(n) = 6, 

V^n{n-6) \ _ 

3 ei? / ’ 


(32) 

(33) 


where (•) is the average defined in Sect. 

One can specify the conditions under which the larger 
system behaves as a reservoir of sides and curvatures for 
the smaller one: suppose that A^2 is the larger system {i.e. 
containing more degrees of freedom) and is the smaller 
one. By definition, the values of effective temperature and 
chemical potential of a reservoir should not be affected by 
the presence of the smaller system, that is: 


dS2 

dK2 

dSo 



fdS2 ' 

\ 


K2 = Kiu2-Ki 

[dK2, 

)k2=kV 

(34) 

N2—Ni\j2—Ni 


N2=N2 



fdS2 " 

\ 


K2=Kiu2-Ki ~ 

\dN 22 

Ik2=kI’ 

(35) 

A^2—-^^1U2— 


N2=N2 



As shown in Appendix C, the derivatives dK 2 S 2 {K 2 , N 2 ) 
and dN 2 S 2 {K 2 , AI 2 ) involve the quantities N 2 /AA 2 and 
K 2 /Af 2 - Therefore conditions (34) and (35) require 


^0 


^0 


N 2 /N 2 — N 2 /N 2 and K 2 /Af 2 ^ K 2 /N 2 when N 2 be¬ 
comes large. Under gentle shuffling, the whole system 
(-5 ^iu 2) remains clustered, so 


N 


1U2 


= 6(A/i -I- A/2 + 1 — x) — c\/Afi -I-A/2 


(36) 


and 

Ki\j2 = —d\/Afi -f A/2. (37) 

On the other hand, the number of configurations where 
the A/i bubbles are dispersed when the two systems are 
mixed is much larger than the number of configurations in 
which remains clustered, as is illustrated in Fig. ([^. 
Thus, A^i ~ A/i and Ki ^ A/). Finally, 5^2 constitutes 
a cluster with A/i holes: N 2 = A^iu2 — A^i, and K 2 = 
Kiu 2 — Ki. Therefore, the condition under which S ^2 acts 
as a reservoir for is: 

\/A/2 Afi , (38) 

what is a more stringent condition than the one (A/2 
Afi) required for additive quantities {i.e.: with no surface 
contributions) |26| . 



Fig. 8. Illustration of two foam samples 5^i and 5^2 S> that 
are mixing together. At the initial stage, (grey bubbles) 
was clustered. Because of the agitation, bubbles do not remain 
clustered but scatter in 5^2 (white bubbles). d /2 is still clustered 
but contains “holes”. 
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6.3 Grand-canonical description 


One can now legitimately use a grand-canonical ensem¬ 
ble to describe statistics of bubble shape. Without loss of 
generality, one considers that is composed of a single 
bubble, with given size i? (alternatively, one could pick 
up a small subset of bubbles within the foam. That is the 
procedure that has been used in [T]). Thus, 




if Ny = n and Ky = Ki.b(i?, n), 
otherwise. 


(39) 


In agreement with Eqs (241 and (251, the probability for 
the specified bubble to have n sides is: 


P(n\R) 




(40) 


0.8 

0.7 


(a) small bubbles 


0.7 

0.6 


(b) medmm bubbles 

|! Potts 



where Ky = K — Krb(7?,n) and Ny = N — n are the 
respective curvature and number of sides of the reservoir 
consisting of the Af — 1 other bubbles. The entropy 
of the reservoir associated with the topological degrees of 
freedom {m} is defined as Sy = lnl7^\ Sy (K^, N^i) 
can be expanded near {K,N). As shown in Appendix C, 
successive terms in the expansion are smaller by a factor 
1/Af, as for additive variables |26]. Therefore: 


Fig. 9. Theory vs Simulations: conditional probability P{n\R) 
for (a) small, (b) medium, and (c) large bubbles of a tridisperse 
foam (size ratios and number fractions are identical to those of 
Fig. (d) distribution of number of sides P{n) in the whole 
foam. 


7 Conclusions and perspectives 


Ssi {Ky = K - Krb, = N - n) = Sy{K, N) 
dN^Sy{K,N)n + 0{M- 


(41) 


and thus, in the limit Af —>■ oo: 

{Ky, Ny) = d§ {K, TV) (42) 


with 


and 


= hm 

A/”—>-oo \dK^ J 


(43) 


u = — lim 


(44) 


\dNet) Kse=P' 

Nse=N 

I3~^ and are the effective “temperature” and “chem¬ 

ical potential” of the reservoir, respectively. 

Eq. (401 eventually leads to the desired expression (Q 
for the conditional probability P{n\R). It can be noticM 
that one recovers the distribution intuited independently 
by Schliecker and Klapp ElEI and Sherrington and cowork¬ 
ers but here with an explicit dependence on the 

bubble size R. This dependence is necessary to reflect the 
correlations between topological and geometrical proper¬ 
ties. 

To illustrate the accuracy of the model. Figure 
compares the distributions predicted by the model with 
those obtained from Potts simulations, for the same tridis¬ 
perse foam as in Fig. The agreement is good for the 
conditional probabilities P{n\Ra) {a G {1,2,3}), and 
especially for the distribution of side number P(n) = 
E^P(R^)P(nlRa). 


In summary, the different assumptions on which the statis¬ 
tical model is based have been systematically stated and 
argued in the present study, providing solids foundations 
to this theoretical framework. 

This study also provides insights on possible improve¬ 
ment and extension of the model: foam coarsening could 
be taken into account to study the evolution of disorders 
on long time scales. The model can also be extended to 
biological tissues; in that case, shuffling can arise natu¬ 
rally from cell activity. Energy must be adapted, and cell 
division and apoptosis must be taken into account. 

Extension to 3D foams can be envisaged, but is not 
straightforward: unlike in 2D foams, the total number of 
sides (nor the total number of faces) of a foam that tiles 
the entire space is not imposed by the space-filling con¬ 
straints in the limit JV ^ oo. Euler’s rule combined with 
Plateau’s laws just allow to relate the mean number of 
faces (/) to the mean number of sides (n) of a bubble [5]: 


(/) 


12 

6 — (n) ’ 


(45) 


(/) and (n) depend on the foam polydispersity [55]. More¬ 
over, a T1 event in a 3D foam is the transformation of a 
triangular face into a Plateau border. Consequently, not 
all the faces of a bubble have the opportunity to be in¬ 
volved in a Tl, and the number of sides or faces are not 
invariants of a 3D shuffled foam. 
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Appendix A: fluctuations 

Fluctuations of bubble areas, pressures, and curvatures 

Strictly, the pressure and area of a bubble fluctuate with 
configurations in a shuffled foam: only the number of gas 
molecules that it contains is conserved (in absence of foam 
coarsening). The variations of area and and pressure as¬ 
sociated with a T1 event are imposed by the space Ailing 
constraints: consider for instance the transition from a reg¬ 
ular, hexagonal foam to a new state obtained through a 
single T1 process. Gauss-Bonnet formula [Eq. require 
that the sides of the 5 and 7 cells are curved, which imply 
that the pressures in the bubbles are not identical any¬ 
more. A T1 event is not an isochoric or isobaric process. 
However, in simulations it is often easier to consider that 
bubble areas are preserved (incompressible gas). It is thus 
implicitly understood that the content of gas in the bub¬ 
bles is not preserved in order to satisfy Young-Laplace and 
Plateau laws. 

Pressure and area of any bubble i are related through 
the equation of state of the gas, so their fluctuations 
have the same order of magnitude: typically SUi/IIi = 
—gSAijAi with g = 0(1) (e.g. g = 1 for an isothermal 
process, g = 3l2 for an adiabatic process when bubbles 
contain a monoatomic ideal gas). However, these fluctu¬ 
ations are generally very small: SAi/Ai ^ SUi/IIi 1, 
and a bubble can be unambiguously identified by its area 
(which is then indistinguishable from its mean value A/). 
The criterion to discriminate the sizes of two bubbles i and 
j is that \Ai — Aj\l(Ai -I- Aj) is much larger than SAi/Ai 
and SAj/Aj. 

It must be emphasized that although the pressure fluc¬ 
tuations are generally negligible (SUi/Ui <C 1), the fluc¬ 
tuations or the pressure difference between two neighbour¬ 
ing bubbles i and j can be as large as | (SUi — SUj)/ (Ui — 
nj)\ ~ 0(1). The variation is specially important when 
one of the bubbles, say bubble i, is directly involved in a 
T1 event. In the frame of our mean field approximation, 

Ui - Uj "fKrh(Ri, rii) = 7^ . (46) 

Hence, the relative variation is ~ 1/ (ni—6), as rii increases 
or decreases by one. 


Fluctuations of K and N 

The fluctuations of the total curvature K and total side 
number iV of a set of JV bubbles depend on the specific 
boundary conditions: for an unbounded foam, both K and 
N are fixed and the fluctuations are trivially zero. 

For a set of Af individual bubbles taken in a larger 
foam, the application of central limit theorem on these Af 


independent random variables yields: SN/N ~ SK/K ~ 

When the Af bubbles are clustered, their side numbers 
an bubble curvatures are independent variables, because 
of the space-filling constraints. However, one can obtain 
the desired fluctuations with the following argument: let 
us “cut” the outer boundary into b (b 1) pieces. For a 
large cluster, the number of sides of each one of them are 
independent random variables, and then, using the central 
limit theorem again, the fluctuations of the total number 

of sides of the boundary W satisfy SN^/Nb ~ W 
Similarly, the fluctuations of the total curvature follow 

SK/K ~ Nb For a round cluster, Nb Af^^'^, and 
thus SK/K ~ The total number of sides of the 

cluster, N, and those belonging to the outer boundary, Nb 
are related by iV = 1!)(N +l — x) — Nb- Therefore, for a large 
cluster (Af ^ 1), SN/N ~ The decays of SK/K 

and SN/N are then different from those corresponding to 
a set of individual bubbles, and correct the scalings given 
in [T]. The scaling SN/N ^ also holds for a foam 

enclosed in a container (while SK/K = 0). In all cases, 
fluctuations become negligible as Af —?> oo. 


Appendix B: Cellular Potts model 

The cellular Potts (or extended large-Q Potts) model has 
been applied to studies in several fields of physics or biol¬ 
ogy, among which: grain growth in metals, or coarsening 
in foams [53], cell sorting m, or foam under shear [36] . 
The model discretizes the continuous cellular pattern onto 
a lattice, with an integer “spin” ak defined at each lattice 
site k, chosen from {!,..., Af}. Thus spin values merely 
act as labels for the Af bubbles. Each bubble extends over 
many lattice sites: typically, each bubble in the simulations 
presented here covers approximately 300 lattice sites. The 
Hamiltonian of the system reads: 

H = + | (47) 

sites bubbles * 

{k,l) i 


The first sum is carried over neighbouring sites (k, 1) and 
represents the surface energy: each pair of neighbours hav¬ 
ing unmatching indices determines a bubble wall and con¬ 
tributes to the bubble wall surface energy. The surface ten¬ 
sion 7 is then proportional to J, the coupling strength be¬ 
tween neighbouring sites, with a prefactor which depends 
on the number of neighbours considered in the coupling. 
In the simulations presented here, the energy is evaluated 
with fourth nearest neighbour interactions (so each lattice 
pixel interacts with 20 other pixels) to smooth the effect 
of lattice anisotropy |28] . 

The second sum in Eq. (47), carried over all the bub¬ 
bles that constitute the foam, is the compressive energy 
of the bubbles. B is the bulk modulus of the gas, and A/ 
is the target area of the bubble i, .i.e the area that the 
enclosed gas would occupy if its pressure were equal to the 
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surrounding pressure. It can be noticed that the expression 
of the elastic term differs slightly from the one commonly 
used in literature mMMi, because of the term in 
the denominator. This makes no difference for monodis- 
perse foams (just a rescaling of the bulk modulus). Here, 
one must take the exact expression of compressive energy 
for bubbles with different sizes. The system evolves using 
Monte Carlo dynamics. The algorithm slightly differs from 
the standard Metropolis algorithm: a site is chosen at ran¬ 
dom, but its value is reassigned only if it is at a bubble 
wall, and then only to one of its unlike neighbours. This 
has two effects: one disregards unrealistic states in which 
bubble walls spontaneously emerge inside a bubble, and 
speed up computations. The algorithm used also forbids 
bubble fragmentation and satisfies the detailed balance 
equation |^. The probability of accepting the trial reas¬ 
signment is P = 1 if the associated change of energy AH 
is negative, and P = exp{—A'H/0) if AH > 0. O defines 
the of the Monte Carlo simulations. 

Parallel tempering technique has been used for 
simulating the foam at different values of temperature si¬ 
multaneously. This technique allows to sample the foam 
dynamics efficiently at very low temperature. 


Appendix C: microcanonical description 


In this Appendix, the microcanonical entropy of a shuf¬ 
fled foam above the cristallization threshold is derived, 
. The expressions obtained from the grand-canonical de¬ 
scription are recovered. Although a little bit trickier, the 
microcanonical description allows to establish the scal¬ 
ing of the successive derivatives of the entropy with Af, 
which is useful to give a quantitative definition of a reser¬ 
voir [Eq. (381] and to justify the expansion used in the 
grand-canonical description [Eq. ©I- As in the Grand- 
canonical description, only needs to be evaluated to 
study the correlations between topological and geometri¬ 
cal features. One treats bubbles of same size as discernable 
objects, the indiscernability factor being incorporated in 
Qis), For simplicity, one assumes that = 3. According 
to self-consistency equations (211 and (22), every accessi¬ 
ble configuration must satisfy: 


At 


= y/Af(P)<5, 


At 


(48) 


with Xi = — 3)/^/3eRj,, S = — 

^Af)/\/3eJ\f(R) and = 3y/TrN'{R~^)/e + K. To get 
an analytic approximation of one treats the Xi (or 
Hi) as continuous variables. Equalities ( |4^ are respec¬ 
tively the equations of a Af—hyperplane lying at a dis¬ 
tance 6 from the origin O, and of a Af—hypersphere of 
radius TZ centered on O. Since K ~ —'/N and N > 
^ 1, one has 0 < S < TZ. The intersection of these 
two surfaces is the {JV — 1)— hypersphere with radius 
p = VTZ'^ — and centered on P, the projection of O 
on the hyperplane (see Fig. [^. Its surface is Sj^-i{p) = 



Fig. 10. Illustration in three dimensions (Af = 3) of the vol¬ 
ume occupied by the accessible states (in red) in the a^i-space. 


27r(A^-i)/2/^-2/r[(Af- l)/2], where P is the Euler’s 
Gamma function. S'_A/_i(p)dpd<5 is the hypervolume of the 
region defined by the intersections of the hypershell of 
thickness dp = TZdJZjp = dAr/(2p) with the two hyper¬ 
planes lying at a distance b and 5-|-d^ from O respectively, 
with d5 = !yJZeAf {R) ■ Within the continuous ap¬ 

proximation, the number of states fZ^*'\K,N\dK,dN) 
whose total curvature and number of sides lie between 
[at. A"-b dAr[ and [Af, A^-|-dAf [, respectively, is simply 
the ratio of this volume to the hypervolume occupied by 
one microstate, v = Jlti with bxi = tt i/VV3i7^- It 
yields: 

n^^\K,N-dK,dN) = 

gAf(in R)/2 (3ev^)^/^p^ dp dS 

P (^^) P P ' 

dp/p and di5/p are slowly varying functions of Af. Using 
Stirling approximation for P when Af S> 1, the expression 
of the entropy of the foam S{K, N) = In 17^*) simplifies to: 


5(i^,iV):=.^|ln 


N 

Af 



-3 (A)-M-b6ee^ -b (InA) W (50) 


where e is the Euler’s number. In the limit N —)■ oo, K/Af 


and N/Af tend to constant values (see Sect. 6.2), and the 


entropy becomes extensive. In particular, for a reservoir 
(that is, a large bubble cluster), K/Af —>■ 0 and N/Af —)■ 6, 
and the entropy reads: 

= y {ln[l8eV^((A-i)-(A)-i)] +(lni?)}. (51) 
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One can now justify the limitation up t o fir st order 
used in the grand-canonical description [Eq. (411] for the 


expansion of the entropy of the reservoir Sag [K^, N^) = 
In 17^ {K^, N^) near (itl, N)-. using Eq. (501, it is easy to 
see that, as A/” —t oo, 


V dK^dN^ J 

Nm=N 




(52) 


for any integers p,q> 0. Therefore, successive terms in the 
expansion are smaller by a factor 1 /M, as it is for additive 
quantities. Note that such a scaling is not mandatory when 
the entropy depends on non-additive quantities as it is the 
case here. In particular, the first derivatives are: 


\ 3ee 

Q 


dKag J Kse = K 
Nse=N 

dS^ \ ^_(E_o 
ONsiJk^^k J {R)Q 

Nse=N 


2e, 


with 


Q = 2ev^^9(i?-i)- 


N 

Jf 


-3 {R)-^ -h6ee 


(53) 

(54) 

K 

... 

The two derivatives become intensive quantities in the 
limit A/" —t 00 . In agreement with Eqs (431 and (441, 
the asymptotic values of the two derivatives define the 
temperature and chemical potential of the reservoir. The 
reader can easily check that their expressions are identi¬ 
cal to those obtained from the grand-canonical description 
[Eqs ^ and (|^] above the crystallization threshold. This 
is not a surprise, as both derivations use the same ap¬ 
proximation, in which the {uj} are treated as continuous 
variables HE]. 

The grand partition function r: is related to S by: 


In = S — jSK + fiN. 


(56) 


Hence, maximization of S under fixed values of K and N 
is formally equivalent to the minimization of the Grand 
Potential J = Ins' (under fixed values of (3 and /i). 

Note that the continuous approximation requires that, 
for every f, 5xi is much smaller than hi and TZ, i.e. Af ^ 1 
and Af{R~^) ^ 7?“^. These conditions are ensured when¬ 
ever there is a large number of bubbles with same size. 
Moreover, the number of accessible states 17^*^ must be 
large, which requires ISeTr ((i?“^) — (i?)“^) > 1 for 

large Af. This is a lower bound to the range of polydis- 
persity on which the continuous approximation is valid. 
Indeed, at lower dispersity, the model predicts a crystal¬ 
lization transition |3]. There is also an upper limit to the 
range of validity of the continuous approximation: most 
of the intersection between the two hypersurfaces must 
lie in the space domain {xi > 0,i = 1,..., A/”} to satisfy 
the constraints Hi > 3. The proportion of the intersection 
lying in this space domain is not easy to calculate. Numer¬ 
ically, it has been established that the condition is satisfied 
whenever {R){R~^) — 1 0.35 and ((T„/(n))^ < 0.4 |3]. 
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